---
title: "R14.Evora_Occupations"
---

# Load Packages
```{r}
library(tidyverse)
library(sandwich)
library(stargazer)
```

# Load Data
```{r}
load("Data/data_final.Rdata")
load("Data/gnr_panel.Rdata")
load("Data/ideo_list.Rdata")
```

# Restrict to Evora
```{r}
data_final_7 <- data_final %>% filter(DistN == 7) %>% filter(!is.na(Evora_TotArea)) 
```

# Run Linear Probability Models Where High Land Reform Saturation (exp_cat2 == 2) is the Outcome (App. Table 4e)
```{r}

models<-list()
models[[1]] <- lm(exp_cat2 == 2 ~ Evora_PcOccExp  + 
                    ln_PopDens_70 + Cities_10km + ln_Confs_All + PrRecH_Adm3 + Isolados1950 + PrWheatAr_1972 + LandIneqRat + TerrainSlopeIndex + MedAltitude + Communist_1975 +  MinDistPCP_Map, 
                  data = data_final_7)
models[[2]] <- lm(exp_cat2 == 2 ~ Evora_POccOnly + 
                    ln_PopDens_70 + Cities_10km + ln_Confs_All + PrRecH_Adm3 + Isolados1950 + PrWheatAr_1972 + LandIneqRat + TerrainSlopeIndex + MedAltitude + Communist_1975 +  MinDistPCP_Map,
                  data = data_final_7)
models[[3]] <- lm(exp_cat2 == 2 ~ Evora_OOnly_pAff + 
                    ln_PopDens_70 + Cities_10km + ln_Confs_All + PrRecH_Adm3 + Isolados1950 + PrWheatAr_1972 + LandIneqRat + TerrainSlopeIndex + MedAltitude + Communist_1975 +  MinDistPCP_Map,
                  data = data_final_7)

# Robust standard errors
robust_se1 <- sqrt(diag(vcovHC(models[[1]], type = "HC2")))
robust_se2 <- sqrt(diag(vcovHC(models[[2]], type = "HC2")))
robust_se3 <- sqrt(diag(vcovHC(models[[3]], type = "HC2")))

# Results
stargazer(models,
          se = list(robust_se1, robust_se2, robust_se3),
          type = "text",
          keep = "Evora",
          covariate.labels = c("% Parish Land Occupied and Expropriated", "% Parish Land Occupied", "% Land Reform Affected Land Occupied"), 
          dep.var.labels = "High Land Reform Saturation", 
          add.lines = list(c("Covariates", "Yes", "Yes", "Yes")))

```